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' Abstract. The Painleve test is very useful to construct not only the Laurent-series 

solutions but also the elliptic and trigonometric ones. Such single-valued functions 
are solutions of some polynomial first order differential equations. To find the elliptic 
solutions we transform an initial nonlinear differential equation in a nonlinear alge- 
braic system in parameters of the Laurent-series solutions of the initial equation. The 
number of unknowns in the obtained nonlinear system does not depend on number 
of arbitrary coefficients of the used first order equation. In this paper we describe the 
corresponding algorithm, which has been realized in REDUCE and Maple. 
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1 INTRODUCTION 



The investigations of the exact special solutions of nonintegrable systems play an important role 
in the study of nonlinear physical phenomena. There are a few methods to construct such solu- 
tions [1-6] in terms of rational, hyperbolic, trigonometric or elliptic functions. These methods use 
the results of the Painleve test, describing behavior of solutions in the neighbourhood of their 
singular points, but do not use local solutions obtained as the Laurent series. In 2003 R. Conte 
£ — ' and M. Musette [7] have proposed the method, which uses such solutions. This method constructs 

global single- valued solutions in two steps. The first step is construction of the local special solutions 
as the Laurent series. The second step is construction of the first order polynomial autonomous 
differential equations which have the same Laurent series solutions. The general solutions of these 
equations are special solutions of the initial system. In this paper we present the computer algebra 
program, which realizes this method. 



2 THE PAINLEVE ANALYSIS 

When we study some mechanical problem the time is assumed to be real, whereas the integrability 
of motion equations is connected with the behavior of their solutions as functions of complex time. 
Solutions of a system of ordinary differential equations (ODE's) are regarded as analytic functions, 
maybe with isolated singular points. A singular point of a solution is said critical (as opposed to 
noncritical) if the solution is multivalued (single-valued) in its neighbourhood and movable if its 
location depends on initial conditions [8]. The general solution of an ODE of order N is the set of 
all solutions mentioned in the existence theorem of Cauchy, i.e. determined by the initial values. 
It depends on N arbitrary independent constants. A special solution is any solution obtained from 
the general solution by giving values to the arbitrary constants. A singular solution is any solution 
which is not special, i.e. which does not belong to the general solution. A system of ODE's has 
the Painleve property if its general solution has no movable critical singular point [9] . 

The Painleve test is any algorithm, which checks some necessary conditions for a differential 
equation to have the Painleve property. The original algorithm, developed by P. Painleve and used 
by him to find all the second order ODE's with Painleve property, is known as the a-method. 
The method of S.V. Kovalevskaya [10] is not as general as the a-method, but much more simple. 
The remarkable property of this test is that it can be checked in a finite number of steps. This 
test can only detect the occurrence of logarithmic and algebraic branch points. To date there 



is no general finite algorithmic method to detect the occurrence of essential singularities 1 . In 
1980, developing the Kovalevskaya method further, M.J. Ablowitz, A. Ramani and H. Segur [12] 
constructed a new algorithm of the Painleve test for ODE's. This algorithm appears very useful to 
find solutions as a formal Laurent series. First of all, it allows to determine the dominant behavior 
of a solution in the neighbourhood of the singular point to- If the solution tends to infinity as 
(t — to) 13 , where (3 is a negative integer number, then substituting the Laurent series expansions 
one can transform nonlinear differential equations into a system of linear algebraic equations on 
coefficients of the Laurent series. All solutions of an autonomous system depend on the parameter 
to, which characterizes the singular point location. If a single- valued solution depends on other 
parameters, then some coefficients of its Laurent series have to be arbitrary and the corresponding 
systems have to have zero determinants. The numbers of such systems (named resonances or 
Kovalevskaya exponents) can be determined due to the Painleve test 2 . 

In [7] the following classical results have been used to construct the suitable form of the first 
order autonomous equation: 

1. The Painleve theorem [9]. Solutions of the equation 

P(y{t),y t (t),t)=0, 

where P is a polynomial in both y(t) and y t (t) = ^jp- has no movable essential singular point. 

2. The Fuchs theorem [14]. If the equation 

m 
fc=0 

where Pk(y(t), t) are polynomials in y{t) and analytic functions in t, has no critical movable singular 
points, then the power of Pk(y) is no more than 2m — 2k, in particular, P m (y) is a constant. 

Therefore, the necessary form of a polynomial autonomous first order ODE with the single- 
valued general solution is 

m 2m — 2k 

2 S a jkV j Vt=^, a 0m = l, (1) 

k=0 j=0 

in which m is a positive integer number and ajk are constants. 

3. The Briot and Bouquet theorem [15]. If the general solution of a polynomial autonomous first 
order ODE is single- valued, then this solution is either an elliptic function, or a rational function 
of e' ,x , 7 being some constant, or a rational function of x. Note that the third case is a degeneracy 
of the second one, which in its turn is a degeneracy of the first one. 



3 THE ALGORITHM AND ITS REALIZATION 

3.1 The Laurent-Series Solutions 

To analyze the method of the Laurent series solutions construction let us consider the generalized 
Hcnon-Heiles system with an additional non-polynomial term, which is described by the Hamilto- 
nian: 

H = g (ft + Vt + ^ 2 + X 2 y 2 ) + x 2 y - - y 3 + JL 
and the corresponding system of the motion equations: 

x tt = -\ lX -2xy+±, 
Vtt = -\2y-x 2 + Cy 2 , 

1 Different variants of the Painleve test are compared in [11, R. Conte paper] 

2 In such a way we obtain solutions only as formal series, but for some nonintegrable systems, for example, 
the generalized Henon-Heiles system [13], the convergence of the Laurent- and psi-series solutions has 
been proved. 
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where x tt = and y tt = Ai, A 2 , [i and C are arbitrary numerical parameters. Note that if 
A 2 7^ 0, then one can put A 2 = sign(\<i) without loss of generality. If C = 1, Ai = 1, A 2 = 1 and 
yu = 0, then (2) is the initial Henon-Heiles system [16]. 

The function y, solution of system (2), satisfies the following fourth-order equation, which does 
not include \i: 

20(7 

y tm = (2C - 8)y tt y - (4Ai + A 2 )y tt + 2{C + l)y 2 t + — y 3 + (4CAi - 6A 2 )y 2 - 4A 1 A 2 y - AH. (3) 

We note that the energy of the system H is not an arbitrary parameter, but a function of initial 
data: y , y ot , y ott and y ottt . The form of this function depends on /i: 

H = \{yl t + Vl) - U + ( X ~k + Vo) (Cyl - A 2% - y ott ) + ^ + ^ - ^ + M 



3- ' V2 ' "V — — — 2(C»g - A 2% - %tt ) • 

This formula is correct only if x = Cyf t — A 2 y — y ott ^ 0. If x = 0, what is possible only at 
yu = 0, then wc can not express x Qt through y , y Qt , y Qtt and y ottt , so H is not a function of the 
initial data. If y ottt = 2Cy y ot — A 2 j/ 0t , then cq. (3) with an arbitrary H corresponds to system (2) 
with fi = 0, in opposite case eq. (3) does not correspond to system (2). 

The Painleve test of eq. (3) gives the following dominant behaviors and resonance structures 
near the singular point t^: 

1. The function y tends to infinity as 6_ 2 (t — to) -2 , where 6_ 2 = —3 or fr_ 2 = ^. 

2. For fr_ 2 = —3 (Case 1) the values of resonances are 



r = -1, 10, (5 ± Vl-24(l + C))/2. 



In Case 2 (6_ 2 = £) 



r = -1, 5, 5± v 7 ! -48/C. 



The resonance r = — 1 corresponds to an arbitrary parameter to- Other values of r determine 
powers of t (their values are r — 2), at which new arbitrary parameters can appear as solutions 
of the linear systems with zero determinant. For integrability of system (2) all values of r have 
to be integer and all systems with zero determinants have to have solutions at any values of free 
parameters included in them. It is possible only in integrable cases. 

For the search for special solutions, it is interesting to consider such values of C, for which r arc 
integer numbers either only in Case 1 or only in Case 2. If there exist negative integer resonances, 
different from r — — 1, then such Laurent series expansion corresponds rather to singular than 
general solution. We demand that all values of r, but one, are nonnegative integer numbers and all 
these values are different. From these conditions we obtain the following values of C: C = — 1 and 

C = -4/3 (Case 1), or C = -16/5, C = -6 and C = -16 (Case 2, a = 1 ~ %/l ~ 48/C ), and also 
C = —2, in which these two Cases coincide. It has been shown in [17] (for /j, = 0) and [6] (for an 
arbitrary value of fi) that single- valued three-parameter special solutions exist in two nonintcgrable 
cases: C = —16/5 and C = —4/3 (Ai and A 2 are arbitrary). 

When the resonance structure is known it is easy to write the computer algebra program, which 
finds the Laurent series solutions with an arbitrary accuracy. For example, we have found 65 coef- 
ficients of the Laurent series for both above-mentioned values of C, the sizes of the corresponding 
output files are about 10 Mb. 



3.2 Two Methods for Construction of Global Single- Valued Solutions 

We have found local single- valued solutions. Of course, the existence of local single- valued solutions 
is a necessary, but not a sufficient condition for the existence of global ones, because solutions, 
which are single- valued in the neighbourhood of one singular point, can be multivalued in the 
neighbourhood of another singular point. So, wc can only assume that global three-parameter 
solutions are single- valued. If we assume this and moreover that these solutions are elliptic functions 
(or some degenerations of them) , then we can seek them as solutions of some polynomial first order 
equations. 
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The classical method to find special analytic solutions for the generalized Henon-Heiles system 
is the following: 

1) Transform system (2) into eq. (3). 

2) Assume that y satisfies some first order equation, substitute this equation in (3) and obtain 
a nonlinear algebraic system. 

3) Solve the obtained system. 

The second way proposed by R. Conte and M. Musette, is the following: 

1) Choose a positive integer m and define the first order ODE (6), which contains unknown 
constants djk- 

2) Compute coefficients of the Laurent series solutions for (2) or (3) with some fixed C. The 
number of coefficients has to be greater than the number of unknowns. 

3) Substituting the obtained coefficients, transform eq. (8) into a linear and overdetermined 
system in djk with coefficients depending on arbitrary parameters. 

4) Eliminate the ajk and obtain the nonlinear system in five parameters. 

5) Solve the obtained system. 

To obtain the explicit form of the elliptic function, which satisfies the known first order ODE, 
one can use the classical method due to Poincare, which has been implemented in Maple [19] as 
the package "algcurves" [18]. 

The second way has a few preferences. The first preference is that one does not need to transform 
system (2) to one differential equation cither in y or in x. Moreover at C = —16/5 not x, but x 2 
may be an elliptic function. To construct the Laurent series for x 2 is easier than to find the fourth 
order equation in x 2 . The main preference of the second method is that the number of unknowns in 
the resulting algebraic system does not depend on number of coefficients of the first order equation. 
For example, eq. (6) with m = 8 includes 60 unknowns djk, and it is not possible to use the first 
way to find similar solutions. Using the second method we obtain (independently of the value of 
m) a nonlinear algebraic system in five variables: Ai, A2, H and two arbitrary coefficients of the 
Laurent-series solutions. 

The first way also has one important preference, ft allows to obtain solutions for an arbitrary 
C, whereas using the second method one has to fix value of C to construct the Laurent series 
solutions, because the resonance structure depends on C. 

3.3 The Computer Algebra Algorithm 

Let us consider computer algebra procedures, which assist to construct the first order equation in 
the form (1) with the given Laurent-series solutions: 

00 

y=Yl c ( fc )* fe ' 

k— — p 

where p is some integer number. We can eliminate from eq. (1) terms more singular than y™ : 

m j<=(m—k)(p+l)/p 

F(Vt,y) = Yl a jk y j y* = 0, a 0m = l. (4) 

k=o j=o 

At singular points j/J™ tends to infinity as t m (p +1 ), so we can present F(y t ,y) as the Laurent series, 
beginning from this term: 

jV m „ 

F(y t ,y)= ]T K s t s (5) 

s= — m(p+l) 

and transform (4) in overdetermined algebraic system: K s = in a^ . We choose N max to be more 
than the number of coefficients a^ . The Maple procedures, which make this transformation are 
presented in Appendix. The trivial variant is the following: Procedure quvar(m,p) calculates the 
number of coefficients a^, procedure equa(a 1 m 1 p,yp2,dyp2) constructs the first order equation 
in the form (4) and procedure equalaur(a,m,p, Nmax) constructs the Laurent series (5). Nmax 
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should be more than quvar(m,p). The first computer algebra realization, which generates only 
terms used in future, has been written in AMP [20] by R. Conte. This algorithm bases on the 
a-method of the Painleve test. Our realization bases on transformations of the Laurent series and 
generates only useful terms as well. There exist Maple [19] and REDUCE [21] realizations of our 
algorithm. The Maple realization {procedure equlaurlist(a, m,p, owe, c)} is presented in Appendix. 

Let us consider how this procedure works. We put o(0, m) = 1, other a(i,j) are unknown. The 
procedure equlaurlist(a, m,p, ove, c) does the following: 

1) Calculates Nmax := quvar(m,p) + ove; 

2) Constructs the list which corresponds to eq. (5): fequlist := equalist(a,m,p); For example, 
if m — 2 and p = 2 we obtain 

fequlist := [[o[0, 0], 0, 0], [o[l, 0], 1, 0], [a[2, 0], 2, 0], [a[3, 0], 3, 0], [a[0, 1], 0, 1], [a[l, 1], 1, 1], [1, 0, 2]]; 

3) To simplify the following procedures puts 

Vfc = -m(p + 1).. - p - 1 : c(k) = 0; 

4) Constructs the list of the Laurent series coefficients of F(y t ,y) (laurlist). Coefficient corre- 
sponding t is constructed due to procedure oneequlaur(c 1 fequlist, k,p) as the sum of the corre- 
sponding coefficients of terms a[i , j]y J y{ . These coefficients are calculated by procedure 
monomlaur(c,mon,k,p), where mon — [a[k,j],i,j]. 

The obtained system is linear in a[i, j] and nonlinear in parameters of the Laurent series. This 
system can be transformed into a nonlinear system in parameters of the Laurent series, so the 
number of unknowns does not depend on m. The resulting nonlinear system can be solved using 
the standard Grobner basis method. 

4 CONCLUSION 

The Painleve test is a very useful tool to find single-valued solution. The corresponding computer 
algebra algorithm has been constructed in Maple and REDUCE. The "naive" algorithm calculates 
many terms to be discarded on the following step. Our algorithm calculates only useful terms. 

The author is grateful to R. Conte and V. F. Edneral for valuable discussions. This work has 
been supported by Russian Federation President's Grants NSh-1685.2003.2 and NSh-1450.2003.2 
and by the grant of the scientific Program "Universities of Russia" 03.02.028. 
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Appendix 



quvar : =proc (m: : integer, p: : integer) 

# 10.10.2003 

# This procedure calculates the number of terms in first order autonomous 

# DDE, which solutions tend to infinity as l/t~p. 

# The maximal degree of the derivative is equal to m. 

local k, j , numterm; 

mimterm:=0; 

for k from to m-1 

do for j from while p*j <= (p+l)*(m-k) 

do numterm : =numterm+l ; 

od; 

od; 

return numterm; 
end; 

equa:=proc(a, m::integer, p::integer, yp2, dyp2) 

# 10.10.2003 

# This procedure constructs the first order autonomous polynomial ODE, which 

# solutions tend to infinity as l/t~p. 

# The maximal degree of the derivative dyp2 is equal to m. 
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local equ, k, j, numterm; 
equ:=0; 

for k from to m 

do for j from while p*j <= (p+l)*(m-k) 
do equ := equ+a[j ,k] *yp2~ j*dyp2~k; 
od 

od; 

return equ; 
end; 

equalaur:=proc(a, m::integer, p::integer, Nmax :: integer) 

# 10.10.2003 

# This procedure expands the first order polynomial autonomous ODE in 

# the Laurent series, including terms from l/t~p to t~Nmax. 

# The maximal degree of the derivative dyp is equal to m. 

local max, equ, k, j ,y,dy,equlist,t; 

equ:=equa(a,m,p,yp,dyp) ; 

y:=0; 

for k from -p to Nmax-p do y : =y+c (k) *t~k od; 
dy:=diff (y,t) ; 
max : =quvar (m , p) +1 ; 
if Nmax > max then max:=Nmax fi; 
for k from to m 
do 

dyp(k) : =convert (taylor (eval(dy**k*t~ ( (p+1) *m) ) ,t ,max) ,polynom) 
od; 

for k from to iquo (m* (p+1) ,p) 
do 

yp(k) :=convert(taylor(eval(y**k*t~((p+l)*m)) ,t, max) ,polynom) 
od; 

equlist: = [] ; 

equ : =expand(eval (equ*t~ (- (p+1) *m) ) ) ; 
for k from 1 to max 
do 

equlist : = [op (equlist) , asubs (t=0 , equ) ] ; 
equ:=dif f (equ,t)/k; 
od; 

return equlist; 
end; 



equalist :=proc(a, m::integer, p::integer) 

# 30.10.2003 

# This procedure constructs the first order autonomous DDE is a list. 

# solutions tend to infinity as l/t~p. 

# The maximal degree of the derivative is equal to m. 

local f equlist, k, j; 
f equlist : = [] ; 
for k from to m 
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do for j from while p*j <= (p+1) * (m-k) 

do f equlist : = [op(f equlist) , [a[j ,k] , j ,k] ] ; 
od 

od; 

return fequlist 
end; 



ydegree:=proc(c,n, j ,p) 

# 1.11.2003 

# This procedure constructs the Laurent series for y~n. 

# Solutions tend to infinity as l/t~p. mon: = [a[i, j] ,i, j] . 

local sumy,k; 
if n=l then return c(j) 
else sumy:=0; 

for k from -p to j+p*n 

do sumy : =sumy+c (k) *ydegree (c ,n-l , j-k,p) ; 
od; 

return sumy; 

fi; 

end; 



dydegree:=proc(c,n, j ,p) 

# 1.11.2003 

# This procedure constructs the Laurent series for dy~n. 

# Solutions tend to infinity as l/t~p. mon:=[a[i, j] ,i, j] . 

local sumdy,k; 

if n=l then return (j+l)*c(j+l) 
else sumdy:=0; 

for k from -(p+1) to j+(p+l)*n 

do sumdy : =sumdy+ (k+1) *c (k+1) *dydegree (c ,n-l , j-k,p) ; 
od; 

return sumdy; 

fi; 

end; 



monomlaur:=proc(c,mon, j : : integer, p: : integer) 

# 1.11.2003 

# This procedure constructs the Laurent series for the monomial 

# a[i, j]*y~i*dy*j . 

# Solutions tend to infinity as l/t~p. mon:=[a[i, j] ,i, j] 

local k,coef ,ydeg,dydeg,sum; 
coef :=op(l,mon) ; 
ydeg:=op(2,mon) ; 
dydeg:=op(3,mon) ; 
if ydeg=0 then 
if dydeg=0 then 
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if j=0 then return coef 
else return 

fi; 

else return coef *dydegree (c , dydeg, j ,p) 
fi; 

else if dydeg=0 then return coef *ydegree(c,ydeg, j ,p) 
else sum:=0; 

for k from -p*ydeg to j+(p+l)*dydeg 

do sum:=sum+ydegree(c, ydeg,k,p)*dydegree(c, dydeg, j-k,p) ; 
od; 

return coef*sum; 

fi; 

fi; 

end; 

oneequlaur :=proc(c, fequlist, j::integer, p::integer) 

# 7.11.2003 

# This procedure constructs the j-th term of the Laurent series of the 

# first order autonomous ODE (the list fequlist) . 

# solutions tend to infinity as l/t~p. 

local equj, k, test; 
equj : =0 ; 

test : =nops (fequlist) ; 

for k from 1 to nops (fequlist) 

do equj :=equj+monomlaur(c,op(k, fequlist) ,j ,p) ; od; 
return equj ; 
end; 

equlaurlist : =proc (a,m: : integer ,p : : integer , ove : : integer , c) 

# 10.11.2003 

# This procedure constructs the Laurent series of the 

# first order autonomous ODE with maximal degree of the derivative is equal to m. 

# Solutions tend to infinity as l/t~p. 

# c(k) are the Laurent series coefficients of the function. 

# The length of the resulting list is quvar (m,p)+ove . 

local k, laurlist, fequlist, Nmax; 

Nmax : =quvar (m , p) +ove ; 

fequlist : =equalist (a,m,p) ; 

for k from -m*(p+l) to-p-1 do c(k):=0 od; 

laurlist : = [] ; 

for k from -m*(p+l) to Nmax-m*(p+l) do 

laurlist : = [op (laurlist) .oneequlaur (c, fequlist ,k,p)] od; 
return laurlist; 
end; 
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